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In this paper we study the bilateral filter proposed by Tomasi and Manduchi, as a spectral domain transform 
defined on a weighted graph. The nodes of this graph represent the pixels in the image and a graph signal defined on 
the nodes represents the intensity values. Edge weights in the graph correspond to the bilateral filter coefficients and 
hence are data adaptive. Spectrum of a graph is defined in terms of the eigenvalues and eigenvectors of the graph 
Laplacian matrix. We use this spectral interpretation to generalize the bilateral filter and propose more flexible and 
^ application specific spectral designs of bilateral-like filters. We show that these spectral filters can be implemented 

with /c-iterative bilateral filtering operations and do not require expensive diagonalization of the Laplacian matrix. 
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^ I. Introduction 

The bilateral filter (BF) proposed by Tomasi and Manduchi [l\ has emerged as a powerful tool for adaptive 
processing of multidimensional data. Bilateral filtering smooths images while preserving edges, by taking the 
weighted average of the nearby pixels. The weights depend on both the spatial distance and photometric distance 
which provides local adaptivity to the given data. The bilateral filter and its variants are widely used in different 
applications such as denoising, edge preserving multi-scale decomposition, detail enhancement or reduction and 
segmentation etc. (2)-||6). Bilateral filtering was developed as an intuitive tool without theoretical justification. Since 
then connections between the BF and other well known filtering frameworks such as anisotropic diffusion, weighted 
least squares, Bayesian methods, kernel regression and non-local means have been explored [7||-[[TT|. 
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The BF is data dependent and hence a non-linear and non-shift invariant filter. So, it does not have a spectral 
interpretation in the traditional frequency domain of images. However we would like to have spectral interpretation 
so that we can modify the global properties of the signal by changing its frequency components. To overcome this 
difficulty, we view the BF as a vertex domain transform on a graph with pixels as vertices, intensity values of each 
node as graph signal and filter coefficients as link weights that capture the similarity between nodes. This graphical 
views of the BF is used in (6), J9|, |T2)-[[T4| but spectral design of filters for images using this graph has not been 



studied. 

We can define spectral filters on these graphs where spectral response is calculated in terms of eigenvectors and 
eigenvalues of the graph Laplacian matrix (15)-|[T8|. This spectral interpretation captures the oscillatory behavior 
of the graph signal fT9| and thus allows us to extend of the concept of frequency to irregular domains. This has led 
to the design of frequency selective filtering operations on graphs similar to that in traditional signal processing. 
These graph spectral filters also have a vertex domain implementation. 

In this paper, we interpret the BF as a 1-hop localized transform on the aforementioned graph. Because the link 
weights in the graph are data adaptive, the problem of structure preserving filtering boils down to low pass spectral 
filtering on that graph. We show that the BF can be characterized by a spectral response corresponding to a linear 
spectral decay. We also calculate the spectral response of iterated BF. We extend this novel insight to build better 
and more general bilateral-like filters using the machinery to design graph based transforms with desired spectral 
response. Our design allows one to choose the spectral response of the filter depending on the application which 
offers more flexibility. We also give a theoretical justification for the design using the framework of regularization 
on graphs. We show that these spectral filters do not require computationally expensive diagonalization of the graph 
Laplacian matrix and then provide an efficient algorithm for implementing these filters using the BF as a building 
block. We examine the performance of the proposed filters in a few applications. 

II. Bilateral Filter as a Graph Based Transform 

Consider an input image to the BF. The value at each position in the output image x out is given by the 
weighted average of the pixels in x^ n . 

*out[j] = ^2 x mH (!) 

• l^i W ij 



The weights Wij depend on both the euclidean and photometric distance between the pixels x^ n [i] and x^ n [j]. Let 
Pi denote the position of the pixel i. The weights are then given by 



\\Pi -JPj\\ 2 \ ( (*in[i\ ~ X w [i]) 2 
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w %3 = exp ( j.exp — ] (2) 



Spatial Gaussian weighting decreases the influence of distant pixel^] and intensity Gaussian weighting decreases 
the influence of pixels with different intensities. Intuition is that only similar nearby pixels should get averaged so 
that blurring of edges is avoided. 

Now, consider an undirected graph G = (V, E) where the nodes V = {1, 2, . . . , n} are the pixels of the input im- 
age and the edges E = {(i, j, ?%)} capture the similarity between two pixels as given by the BF weights (Figure[T]). 
Image can be considered as a signal defined on this graph : V — > K where the signal value at each node equals 
the corresponding pixel intensity. 

Adjacency matrix W of this graph is given by W = [wij] nXn . Let D be the diagonal degree 
matrix where each diagonal element Djj = J2i w ij- With this notation the filtering operation ^^X?',/^ 
in ([T]) can be written as (9| ^/^^\ 

X^D^WXin (3) Fig.l:TheBF 

graph 

It can be seen from ([!]) that the output at each node in the graph depends only on the nodes in its 1-hop neighborhood. 
So, the BF is a 1-hop localized graph based transform. Also, note that the BF includes the current pixel in the 
weighted average. So, the graph corresponding to the BF has a self loop i.e. an edge connecting each node to itself 
with weight 1. Other filtering techniques such as Gaussian smoothing and non-local means can also be described 
using similar graphical models J9|, fT2| . 

A. Graph Spectrum and Data Adaptivity of the BF 

Spectrum of a graph is defined in terms of the eigenvalues and eigenvectors of its Laplacian matrix. The 
combinatorial Laplacian matrix for the graph G is defined L = D — W. We use the normalized form of the 
Laplacian matrix given as C = D -1 / 2 LD -1 / 2 . C is a non-negative definite matrix ]20[ . As a result C has an 
orthogonal set of eigenvectors U = {ui, . . . , 112} with corresponding eigenvalues <j{G) = {Ai, . . . , X n }. So, C can 
be diagonalized as 

C = UAU* (4) 

where A = diag{\i, . . . , X n }. 

Similar to classical Fourier transform, the eigenvectors and eigenvalues of the Laplacian matrix C provide a 
spectral interpretation of the graph signals. The eigenvalues {Ai,...,A n } can be treated as graph frequencies, 
and are always situated in the interval [0, 2] on the real line. The eigenvectors of the Laplacian matrix demonstrate 



increasing oscillatory behavior as the magnitude of the graph frequency increases 1 19 1. The Graph Fourier Transform 



! Eq. {1} shows that the averaging is done over all pixels. However in practice, one assumes non-zero weights only for the pixels which 
have \\pi — Pj\\ < 2a s |2(. 
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— graph with bilateral weights 

— graph with gaussian filtering weights 
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Fig. 2: Ek is the fraction of total energy captured by the first k spectral components. 



(GFT) of a signal x is defined as its projection onto the eigenvectors of the graph, i.e., x(Xi) = (x, u^), or in 
matrix form x = U*x. The inverse GFT is given by x = Ux. 

Figure [2] shows the fraction of total signal energy captured by the first k spectral components of a graph signal 
corresponding to a 64 x 64 image block. In this example we consider the following underlying graphs (1) the 
BF graph (2) the graph corresponding to Gaussian smoothing where the weights depend only on the geometric 
distance between two pixels. We can see that due to the data adaptivity of the BF graph, most of the signal energy 
is captured in low frequency part of the BF graph spectrum in comparison with spectrum of Gaussian smoothing 
graph. So, the spectral basis of the BF graph offers better energy compaction for the given signal. 

III. Spectral interpretation of the bilateral filter 
Similar to conventional signal processing, graph spectral filtering is defined as 

x ou t{\) = h(\i)x in (Xi) (5) 

h(Xi) is the spectral response of the filter according to which spectral components of an input signal are modulated. 
Using the definition of GFT and the diagonalized form of £, we can write graph spectral filtering in matrix notation 
as (T7J 

*out = h(A) USq^ = h(C)x in (6) 

^fse Spectral GFT 
response 

To exploit this framework of graph spectral filtering we rewrite the BF in ([3]) as 

X out = D-V2 D -V2WD-V2 D V2x; n 

^B 1/2 x out = (I - £)D 1 / 2 x m (7) 

From this equation, we can see that the BF is a graph transform, similar to the one in ([6]), operating on the 
normalized input signal 5q n = D 1 / 2 x^ n producing the normalized output -k out = D 1 / 2 x ou ^. This normalization 
allows us to define the BF in terms of the non-negative definite matrix C and thus have a spectral interpretation. It 
also ensures that a constant signal when normalized, is an eigenvector of C associated with zero eigenvalue pl| . 
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Fig. 3: Spectral responses of the BF and iterated BF. 



Following ([6]) we have, 



± out = U(I - A)U'x. 



(8) 



This shows that the BF is a frequency selective graph transform with a spectral response fiBF(^i) = 1 — \ which 
corresponds to linear decay (See Figure [3}. The BF tries to preserve the low frequency components and attenuate 
the high frequency components. 

The BF is used iteratively in many applications. There are two ways to iterate the BF (1) by changing the weights 
at each iteration using the result of previous iteration (2) by using fixed weights at each iteration as calculated 
from the initial image. In the first method the BF graph changes in every iteration. So we cannot have a spectral 
interpretation. In the second method the graph remains fixed at every iteration. Here we consider the second method. 
Iterating preserves strong edges while removing weaker details. This type of effect is desirable for applications 
such as stylization EL The BF iterated fc-times can be written in matrix notation as 



Xout = (D^W)* X m = (I - Crf Xj. 



(9) 



where C r = D 1 L is called the random walk Laplacian matrix. It can be shown that any graph transform h(C r 



can be written in terms of C as h(C r ) = D 1 / 2 h(C)'D 1 / 2 |21 Proposition 2]. Using this fact, we can rewrite ([9]) 
as 



Z out = U(I - A)*U*x 



(10) 



The spectral responses corresponding k = 2, 3, 4 are shown in Figure [3] The figure suggests that iterative application 
of the BF suppresses more of the high frequency component which is consistent with the observation. Equations 
([8]) and ( |T0] ) give a different perspective to look at the bilateral filter. They hint at filter designs with better spectral 
responses that can be tailored to particular applications. 



IV. Application specific spectral designs 

The BF and iterated BF have fixed spectral responses. But these responses may not be suitable for all applications. 
Below we discuss two applications (1) image denoising and (2) segmentation to illustrate design of more flexible 
spectral filters. 
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Denoising.: We consider the problem of image denoising with additive zero-mean white noise. 

y[i]=x[i]+e[i] (11) 

where x is the original image that we want to estimate, y is the observed noisy image and e is zero-mean white 
noise with variance a 2 . As explained before most of the signal energy lies in the low frequency part of the BF 
graph spectrum due to its data adaptivity. So most energy in the high frequency spectrum corresponds to noise. We 
can put this intuition in a more principled framework of regularization where the problem of denoising is equivalent 
to minimization of a penalty functional (7J. This penalty functional is composed of two terms. The first term is a 
fit measure and the second term is a data dependent smoothness constraint as captured by regularization operators 



on the BF graph [22], |23|. 



£7(x) = i||y-xf + %\\h p (C)xf (12) 
Note that we normalize x and y as in ([7]). The regularization functional can be written in spectral domain as 

n 

||/i p (£)x|| 2 = ^[/i p (A,)] 2 [£(A,)] 2 (13) 
i=i 

hp(X) > is chosen to be a non-decreasing function in A so that the high frequency components are penalized 
more strongly. Putting dC/d~k = 0we get the optimal x as 

Z opt = (I + phl(C))- 1 y 

= U(I + p^(A))" 1 U t y (14) 

So, for a chosen regularization functional h p (X), the spectral response of the denoising filter is given by 

^ (A) = l + phl(X) (15) 

Denoising filters suggested by the regularization framework are essentially low pass filters in the spectral domain 
of the graph (Figure [5]). The BF has a relatively poor spectral decay profile due to which it blurs the textures in 
the denoising process. 

Image Segmentation.: Shi and Malik [24] formulated the problem of image segmentation as a graph partitioning 
problem on a graph similar to the BF graph. To obtain an m-way partition of the graph, we only need to find the 
projection of the graph signal on the first few eigenvectors with the smallest eigenvalues. This projection gives a 
very coarse version of the signal which then can be used to perform image segmentation. So a suitable filter for 
this application should be a low pass filter in the graph spectral domain with a small cut off frequency and sharp 
transition band (Figure |6j). Iterative application of the BF also gives a coarse version of the image with details 
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removed. But it favours the eigenvector with the smallest eigenvalue |6| which corresponds to the DC component 
of the image as seen from its spectral response. On the other hand, for graph partitioning we need the eigenvector 
with second smallest eigenvalue (and eigenvectors with larger eigenvalues for finer partitions). 

To summarize, different applications require filters having different spectral responses h(X). So, we would like 
to design more general bilateral-like filters with desired spectral response. These filters are of the form h(C) = 
U/i(A)U*. A direct implementation of these filters requires diagonalization of C which is of the order 0(N 3 ). 
For large graphs such as the one considered here, this is computationally very expensive. Fortunately, we can 
approximate any spectral response h(X) by a polynomial in A. These polynomial spectral filters are fc-hop localized 
on the graph where k is the degree of the polynomial. This leads to an easy and efficient implementation scheme 
for these filters as explained in the next section. 



V. Polynomial approximation and fast implementation 

The spectral response of the iterative bilateral filter (Figure |4ja)) given in ( fTO] ) is a degree k polynomial. This is 
a special case of a general class of real polynomials of degree k given as 

k 



h{A)=roH(I-r i A) i 



(16) 



where the roots T{ can be either real or complex conjugate pairs. This generalization allows k + 1 degrees of freedom 
in choosing the spectral response of the filter. Further, the corresponding transform H in pixel domain is a matrix 



polynomial of C r with the same roots as in (16) i.e., 



H = D~ 1 / 2 U/i(A)U t D 1 



/2 



i=l 



This leads to following result: 



(17) 
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(a) (b) 
Fig. 4: (a) Iterated BF (b) Fast implementation of a bilateral-like filter with polynomial spectral response using BF as a building block 



Theorem V.l. Any graph filter on an image having a polynomial spectral response of degree k can be implemented 
in the pixel domain as an iterative k step bilateral filter operation. 



Proof: For k = 1, the filter H = I — r\C = (1 — ri)I + r\D 1 W. The output in this case is given in (18) 
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whose convergence rate is r\ times the convergence rate of the original bilateral filter D _1 W. 

y = Hx = (1 - n)x + nD _1 Wx, (18) 

For k > 1, the transform H is a cascaded form of k such bilateral filtering operations (See Figure |4jb)). ■ 
Further, if h(X) is not a polynomial of A, it can be approximated with a polynomial kernel which can then 
be implemented as a generalized iterative bilateral filtering operation. It has been shown in (15) that minimax 
polynomial approximation of any kernel h(X) not only minimizes the Chebychev norm (worst-case norm) of the 
error between kernel and its approximation, it also minimizes the upper-bound on the error ||i7 exact — i7 a PP rox || 
between exact and approximated filters. In our experiments, we approximate any non-polynomial h(X) with the 
truncated Chebychev polynomials (which are a good approximation of minimax polynomials). 

VI. Examples 

We examine the performance of the proposed spectral design of bilateral-like filters in two applications. First, 
we consider the image denoising problem. We experiment with a spectral response obtained by the regularization 
framework in Section 4. We take the regularization functional h p (X) = A 2 which suggests a denoising filter with 
h(X) = 1/(1 + A 2 ). We take its 5 degree polynomial approximation. Figure [5] shows the denoising results using 
this filter and the BF^j It can be seen that the BF preserves edges, but it blurs the texture in the denoising process 
while the proposed denoising filter does a better job of preserving texture. This is also reflected in the SNR values. 

Next, we consider an iterative application of the BF. Iterated BF removes minor details from the image while 
preserving prominent edges. This can be used as an effective preprocessing step in edge-detection and segmentation 
etc. As stated before, iterated BF favours the DC component which is not useful for segmentation. We use a low 
pass spectral kernel with small cut-off and sharp transition band so that the second (and a few higher) spectral 
components get at least as much weight as the DC component. We use a 20 degree polynomial approximation of 
this kernel and compare it with 20 iterations of the BF. Figure [6] shows that weak edges are blurred more with 
iterated BF (as expected from the spectral response) compared to the proposed filter. 

VII. Conclusion 

In this paper we explained the bilateral filter as a graph spectral filtering operation. With this novel perspective, 
we proposed a family of more flexible bilateral-like filters with desired spectral responses. We gave an easy 
implementation scheme for these filters. Their utility was motivated through few examples. An immediate interesting 
extension to this work would be to explore different spectral filters suitable for particular applications. Another 
topic of interest is the design of filter banks using these bilateral-like filters. 

2 The BF is not the best denoising filter available. We use it in our comparison to emphasize the qualitative differences in filtering results. 




(e) (f) 

Fig. 5: (a) Original image (b) Noisy image, SNR = 20 dB (c) Output of the BF(a r = 0.035, a d = 2), SNR = 20.65 dB (d) Output of the 
proposed filter, SNR = 22.64 dB (e) Spectral response of the BF (f) Spectral response of proposed filter 
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